Field
Embodiments described herein relate generally to a method of reconstructing computed tomography images, and more specifically to reconstructing computed tomography images using more than one data set acquired using more than one detector configuration.
Description of the Related Art
The past several years have seen technological advances that are fueling improvements to X-ray computed tomography (CT), including: new iterative image reconstruction algorithms with novel regularization methods, and rapid development of spectral CT. Specifically, there has been a growing interest in iterative reconstruction algorithms for CT using total-variation (TV) regularization, which can lead to high quality reconstructions from highly undersampled projection data. Simultaneously, the development of spectral CT exhibits numerous advantageous over conventional CT. For example, spectral CT offers the additional clinical information inherent in the full spectrum of an X-ray beam. Moreover, spectral CT facilitates, in discriminating tissues, differentiating between materials such as tissues containing calcium and iodine, and enhancing the detection of smaller vessels. Among other advantages, spectral CT is also expected to reduce beam-hardening artifacts.
One well known iterative technique is the Algebraic Reconstruction Technique (ART). This technique is essentially a method for iteratively solving the matrix equationf=Au, where f is a vector of the projection measurements, u is a vector of the image values, and A is a system matrix corresponding to the discretized Radon transform of the X-ray beams passing through an image object. By recognizing that each row vector {right arrow over (A)}i of the matrix A together with the corresponding projection value fi defines an affine space, an image of the image object can be found through successive affine projections onto the successive affine spaces corresponding to the rows of A. This iterative process converges by using the previous estimate of the image vector um−1i to solve for the current image vector estimate umi using the expressionumi=um−1i−{right arrow over (A)}m−1(fm−1−{right arrow over (A)}m−1·um−1i/{right arrow over (A)}m−1·{right arrow over (A)}m−1),where each iteration progressively estimates umi for m=2, . . . , NData, u10 is the initial guess, and the superscript i indicates the ith iteration of affine projections for all values of m. The iterative process continues until the image estimates converge according to some predefined metric.
Typically, following a series of affine projections, a constraint is imposed in order to ensure that the image u converges to a physically meaningful image. For example, in absorption imaging, the image value must be non-negative because a negative absorption value implies gain, which is not physically realistic. Therefore, the final value after each iteration, uNDatai, is subject to a predefined constraint based on a priori knowledge of the image (e.g., no gain), and the constrained final value is then used as the initial value, u1i+1, for the next iteration of affine projections. Periodically subjecting the image estimates to a predefined constraint is referred to as regularization.
Iterative reconstruction algorithms augmented with regularization can produce high-quality reconstructions using only a few views and even in the presence of significant noise. For few-view, limited-angle, and noisy projection scenarios, the application of regularization operators between reconstruction iterations seeks to tune the final result to some a priori model. For example, enforcing positivity, as discussed above, is a simple but common regularization scheme. Minimizing the “total variation” (TV) in conjunction with projection on convex sets (POCS) is also a very popular regularization scheme. The TV-minimization algorithm assumes that the image is predominantly uniform over large regions with sharp transitions at the boundaries of the uniform regions. When the a priori model corresponds well to the image object, these regularized iterative reconstruction algorithms can produce impressive images even though the reconstruction problem is significantly underdetermined (e.g., few view scenarios), missing projection angles, or noisy.
In addition to the ART iterative image reconstruction method, the reconstruction problem can be solved using the TV semi-norm regularization by using the primal-dual algorithm developed by Chambolle and Pock. The Chambolle-Pock (CP) algorithm applied to the TV semi-norm solves the problem
            min      u        ⁢          {                                    1            2                    ⁢                                                                  Au                -                f                                                    2            2                          +                  λ          ⁢                                                                  (                                                                        ∇                    u                                                                    )                                                    1                              }        ,where the last term, the l1-norm of the gradient-magnitude image, is the isotropic TV semi-norm. This expression is referred to as the primal minimization. The spatial-vector image ∇u represents a discrete approximation to the image gradient. The expression |∇u| is the gradient-magnitude image, an image array whose pixel values are the gradient magnitude at the pixel location. Through a process of first finding the convex conjugate of the expressions inside the minimization function, one can obtain the dual maximization. There are several methods for solving the primal-dual saddle point in order to solve the reconstruction problem. The CP algorithm uses the forward-backward proximal-splitting method to obtain an iterative algorithm that converges to an image satisfying the TV minimization problem. For the case of the TV semi-norm, the pseudo-code for the resulting iterative algorithm is given in Table 1.
TABLE 1Pseudocode for N-steps of the CP algorithm using the TV semi-norm.1:L ← ||(a, ∇)|| 2; τ ← 1/L; σ ← 1/L; θ ← 1; n ← 02:initialize u0, p0, and q0 to zero values3:ū0 ← u04: repeat5:pn+1 ← (pn + σ(Aūn − f))/(1 + σ)6:qn+1 ← λ(qn + σ∇ūn)/max(λ1I, |qn + σ∇ūn|)7:un+1 ← un − τATpn+1 + τDivqn+18:ūn+1 ← un+1 + θ(un−1 − un)9:n ← n + 110: until n ≧ NThe constant L is the    2-norm of the matrix (A, ∇); τ and σ are nono-negative CP alogorithm paramters, which can be set to 1/L; θ ε [0,1] is another CP algorithm parameter, which can be set to 1; n is the iteration index; u0 can be initialized to an initial guess of the image; p0 can be initialized to the projection data; q0 can be initialized to the image gradient; and 1I is a diagonal matrix with ones along the diagonal; and Div is the divergence operator.
Also, rapid advances in energy-discriminating detectors have generated a surge of interest in spectral CT, where the broad x-ray tube spectrum is sampled into spectrally unique channels of data, which could lead to the elimination of common artifacts, patient dose reduction, and new applications of CT.
Here, we adopt notation conventions for two different kinds of matrix norms in order to avoid confusion. The first type is the “entry-wise” p-norm, which for matrix X is defined as
                  X              p    =                    (                              ∑                          i              ,              j                                ⁢                                                                  X                ij                                                    p                          )                    1        /        p              .  For p=2 we get the common “Frobenius norm.”
The second type is the “Schatten p-norm,” which arises when applying the lp-norm to the vector of singular values associated with X. Here the uncommon notation ∥·∥Sp is used, where specifically for a matrix X, the Schatten p-norm is given by
                            X                    Sp        =                  (                              ∑            i                    ⁢                                                                  σ                i                                                    p                          )                    1        /        p              ,where σi denotes the ith singular value of X. Some common special cases are the “nuclear norm” (p=1), the Frobenius norm (p=2), and the “spectral norm” (p=∞).